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Abstract 

The calculation of self-energy corrections to the electron bands of a metal requires the evaluation 
of the intraband contribution to the polarizability in the small-q limit. When neglected, as in 
standard GW codes for semiconductors and insulators, a spurious gap opens at the Fermi energy. 
Systematic methods to include intraband contributions to the polarizability exist, but require a 
computationally intensive Fermi-surface integration. We propose a numerically cheap and stable 
method, based on a fit of the power expansion of the polarizability in the small-q region. We test 
it on the homogeneous electron gas and on real metals such as sodium and aluminum. 
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I. INTRODUCTION 



More than 20 years of successful app 
and its numerical implementations 



)lications have established Hedin's GW approach 
j, |4J as the state-of-the-art and most widely used 
theoretical method for ab-initio bandstructure calculations including self-energy effects. Ef- 
ficient algorithms have been devised to encompass the major numerical bottlenecks in such 
calculations, e.g. by avoiding k-space convolutions by a space-time method [5|, or avoid- 



ing summations over empty states in the evaluation of the po 
localized basis functions [8|, and/or model screening functions 



arizability 

rr 



a, LZJ] , or using 



llOl . Illj . Several computer 

codes have been devised for ab-initio GW calculations, and are presently available, often 



under public domain 



13 



141 ] . However, systems with metallic screening present an addi- 



tional, numerically challenging, difficulty: in the evaluation of the k-space integrals for the 
intraband contribution to the electron screening, the contribution of the Fermi surface can 
dramatically slow down the convergence with respect to the k-space sampling. As a conse- 
quence, the possibility to perform such ab-initio GW calculations in gapless systems with 
a large unit cell is hindered. Even worse, when calculations are performed with standard 
computer codes, unconverged k-point sampling gives rise to a spurious gap at the Fermi 
level. The gap vanishes only in the limit of infinitely dense sampling, and is shown to close 
very slowly as the number of k-points increases. Solutions based on explicit Fermi surface 
integration [l5, 



16 



171 ] are effective but result in cumbersome coding and substantial increase 



of computation time. 

In this paper we present a numerically stable and efficient method, based on a Taylor 
expansion of the polarizability matrix in the small-q region, which includes intraband contri- 
butions and avoids explicit Fermi-surface calculations. The method has been implemented 



successfully into the abinit 12|, ll8j package, and is shown to remove the spurious gap at 



the Fermi level already with a limited number of k-points. Results are presented for the 
homogeneous electron gas (HEG) as well as for real metals such as Na and Al. 

This paper is organized as follows: in Sec. [IT] we briefly review the standard GW scheme 
and describe the difficulties that arise when it is applied to metals naively. In Sec. [TTT] we 
analyze the origin of the problem, and propose our solution in Sec. |TVJ In Sec. [V] we test 
the method on different metallic systems, and we discuss the results in Sec. I VI I 
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II. THEORETICAL FRAMEWORK 



The present work deals with the many-body prob 
based on the following set of self-consistent equations 



em in the standard Hedin's scheme 

I 



G(l,2) = G {l,2)+jG (l,3)E{3,4)G(A,2)d3dA 



(1) 



r(l,2;3) = 5(1,2)5(1,3)+ / S G{4, 6)G{7, 5)r(6, 7; 3) dAd5d6d7 (2) 



X(l,2) 



5G(4,5) 

-i / G(1,3)G(4, l)r(3,4;2) d3d4 



W(l,2) = v c (l,2) + y^(l,3)x(3,4)iy(4,2)d3d4 
£(1,2) = % JG(l,3)W(4 : ,l)r(3,2;A)d3d4, 



(3) 
(4) 
(5) 



G and Go are the exact and Hartree's Green functions for the electron, vc is the bare 
Coulomb interaction, W is the screened potential, x is t ne electric polarizability, E is the 
self energy, and T is the vertex function. An argument such as "1" stands for the set of 
position, time and spin variables (r^, ti, o"i). Equations ([U15]) constitute a formally closed set 
of equations for the five correlators. The functional derivative in Eq. (j2j) provides the vertex 
corrections and is a formidable computational difficulty The most important approximation 
that is usually made is to neglect the vertex entirely and put T(l,2;3) = 5(1,2)5(1,3) in 
the remaining four equations. This explains the name GW, since Eq. (j5J) now simplifies to 
the product 

S(l,2)=zG(l, 2)1^(2,1). (6) 

A GW calculation proceeds as follows. One assumes initially £ = 0, G = Gq in Eq. (JT]). 
Next, one determines x = Xo through Eq. ([3D with T = 1, and computes W from Eq. (J4J). 
The first estimate of £ is obtained in Eq. (JBJ), and can be used to update G and the other cor- 
relators. This procedure can be iterated until self-consistency is reached. However, several 



non-selfconsistent GW approaches are possible [19|. In the present work, we choose to per- 
form calculations within the so-called GqWq approximation 3j, |4j, which stops the iteration 
without updating G and W. One first evaluates the independent-particle polarizability 



Xo(l,2) = -«G (l,2) G (2,l) 



(7) 



3 



and the dielectric function 



e(l,2) 



u c (l,3)xo(3,2)d3d4, 



(8) 



which provides the solution of Eq. (TjJ, Wo(l> 2) = J e -1 (l, 3) vc(3, 2) <i3, and the self energy 
£(1,2) = iG o(l, 2)Wo(2, 1). This approximation is usually a successful one, while self- 
consistent GW has been shown to lead to a worse treatment of electron correlations in 



prototypical systems such as the HEG (where it gives a bandwidth larger than the D 



20| |) and solid silicon (where the band gap turns out even larger than experiment 21]). The 



T one 



simpler GoWo approach leads typically to a 10% bandwidth reduction with respect to DFT, 
in better agreement with experiment, thus suggesting a partial cancellation of errors due to 
lack of self-consistency and of vertex corrections in 22j, [23j . 

In the GqWq approximation, one can start from a DFT-LDA electronic-structure cal- 
culation. The quasiparticle energies Ej are hence evaluated as first order corrections to 
the Kohn-Sham (KS) eigenvalues €j, with respect to the perturbation (£ — Vxc), and by 

,Q: 

(Efo) - V xc ) 



linearizing the energy dependence of £ 



E ~ . 

^3 — 



1 - 



&W \ul=€i 



(9) 



Expectation values are taken on the corresponding KS state |k, j); the denominator is the 
quasiparticle weight. One of the heaviest parts of the GoWo computation is the inversion of 
the symmetrized dielectric matrix, which in reciprocal space reads: 

1 1 



CG,G'(q»w) = 5 g ,g> ~ 4vr- 



X0G,G'iq,^J 



(10) 



|q+G| ' |q+G'| 

The inversion must be performed on a mesh of frequencies spanning a range significantly 
wider than the range of interest for the bandstructure. 

The inverse dielectric matrix leads immediately to the effective screened potential: 

1 n , . 1 



W, 



G,G'iq,wj 



47T- 



iq,^j 



(ii) 



| q+G |-G,G'v-«-/| q+G ,| 

A great simplification can be achieved by introducing an additional plasmon-pole approx- 
imation, where the frequency dependence of each G,G' matrix element is parameterized 
by: 

fi G,c(q) 



; G,G' 



iq,^j 



;g,c 



+ 



^G,G'(q) 



(12) 



The parameters f2 G G , (q) and u) G G / (q) are determined by evaluating the polarizability 
XoG.G't^j u ) om y a t t wo values of the frequency, usually at u = and at a purely imaginary 
frequency of the magnitude of the plasma frequency u = iuip. In the following we adopt this 
plasmon-pole model, since the difficulties related to the small wave-vector screening would 
occur identically if the detailed uj dependence of e -1 were considered. 



The polarizability xo is given by the standard expression 24j, |25] 

2 V- f j3 , /(9'(k + q))-/( e ,(k)) 



X0G,G'(q,^) = -77— I d ^ k 

vbz 7-7 Jbz 

3,3 



iuj - [ej'(k + q) - £j(k)] 
(k,j|e-^+ G > f |k + q,iO(k + q ! j , |e^ +G '> f |k,j) ! (13) 

where /(e) are Fermi occupation numbers at a small smearing temperature, |k,j) are the 
KS states, and the factor 2 accounts for spin. Complex conjugation gives G ,(q, zc<j) = 
Xog'.gCQ; hence also cg,g' is a Hermitian matrix for purely imaginary frequencies. For 
uj 7^ 0, q = and G or G' equal to 0, this expression vanishes exactly because of orthogo- 
nality (j j' terms) or equality of Fermi numbers (j = j'). The rate at which xo vanishes 
as q — > is relevant for contrasting the Coulomb singularity that appears in the dielectric 
matrix. To take care of interband terms (j ^ j'), a standard solution is to expand the matrix 
elements of Eq. (1131) by means of the formula j^g]: 

(Me-^lk + q,/) - < k ^1+^- V -l k '^ + ^/1[^^- g ]l k ^ ,/ ) ; (14 ) 
q^o e j'(k) — e j(k) 

where Vnl is the non local part of the pseudopotential. By substituting this expansion into 
Eq. ( Tl3l) one gets a small q expansion of the polarizability, which can be used to evaluate 
the q — ► limit of ?~ 2 Xo(q) appearing in Eq. ( flOi) when G = G' = 0. Intraband terms 
(j = j') are P u ^ to zero. While this method is satisfactory for semiconductors, it gives rise to 
substantial difficulties for metals, where intraband terms are also important. This leads to 
an incorrect evaluation of the Xo contributions in the q-space region closest to the origin, i.e. 
at one out of N kpt points of the mesh of q points. At first sight, as Xo enters the calculation 
of E through a iVtpt-discretized q-convolution in reciprocal space, one might think that this 
single incorrect value should affect the energy corrections (E — Vxc), with an error of order 
N^.. However, the singular behavior of the Coulomb repulsion vq near q = requires an 

explicit integration around the singular point, which makes the final outcome sensitive to 

_ i 

the incorrect Xoo,o(0,iw) with an error of order A^ k f . 
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Ikl(a ) 



Figure 1: (Color online) GoWo band structure of Na (110 direction), showing the appearance of 
an unphysical gap, and its dependence on different numerical convergence parameters. Panel (a) 
shows the dependence with respect to the number of empty states in Eq. f)13j> : (b) with respect to 
the smearing temperature; (c) with respect to the k-point mesh. Panel (d) shows the dependence 
of the unphysical gap on the inverse number of k-points in each direction; the dashed line is a 
fitted ^N-^ + aiN^f. 

III. NAIVE APPLICATION OF A STANDARD G W CODE TO METALLIC SYS- 
TEMS 

The incorrect small-q values of xo induce the opening of an unphysical gap at the Fermi 
energy in the GqWq band dispersion of simple metals (such as the HEG and sodium), as 
shown in Fig. [IJ The figure also shows the convergence properties of the width of this 
unphysical gap, computed by extrapolation from the two sides. The only significant depen- 
dency is on the number A^ p t of sample points in the k-space mesh: Fig. [T] shows that the 

unphysical gap does tend to close for increasing mesh size, but only extremely slowly, as 
_i 

N kp l , for the reasons discussed at the end of Sec. [TlJ Therefore, it is practically impossible 
to close the gap by brute-force mesh refinement, especially because the computation time 
of the dielectric matrix grows as N£ ' t . 

The spurious gap is essentially independent of most numerical convergence parameters, 
such as the number of empty states and the smearing temperature, as shown in Fig.[T^a,b). A 
larger smearing temperature for electronic occupancy would reduce this unphysical gap, but 
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o HEGr s =3.5 



s < < ■ ' 



b-2 - 



■ ■ DFT-LDA 

■ ■ G W Q band with standard implementation 



Ikl (a ) 



Figure 2: (Color online) Bandstructure for the HEG (r s = 3.5 Oq) computed with the standard 
implementation of the GqWq method (dashed line). The spurious gap, caused by the lack of the 
intraband term in the screening, is removed when the computed polarizability is replaced by the 
Lindhard one (dot-dashed line). The KS band is also displayed for reference (dotted line). 




Figure 3: (Color online) Numerically computed HEG screening function e _1 (q, uj) (r s = 3.5 oo), for 
uj = and uj = iujp, compared to the Lindhard function. For q — > the discontinuously incorrect 
points -pointed at by arrows- appear, due to the lack of the intraband term. The differences at 
large q are due to the finite number of empty states included in the sums of Eq. (p~3|) . 

it is a mere technical device, and convergence should be checked in the limit of vanishingly 
small smearing, where the actual metallic state is recovered. 

The origin of the unphysical gap is the incorrect q = value of the screening function as 
demonstrated in Fig. [5J where the gap is shown to disappear when the numerical dielectric 
matrix is replaced by the Lindhard function [2^]. In metals, the dielectric function e is 
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expected to diverge when both u — > and q — > (by contrast, it goes to its finite static 
limit in semiconductors and insulators). For example for the HEG, interband transitions do 
not contribute to the sum in Eq. (I13p . At the same time, expression (j!4p cannot yield correct 
intraband (j = f) contributions for q = 0, and in practice standard codes do not evaluate 
such terms due to the equality of the occupancy factors. The resulting incorrect null value 
°f Xo(q - ► 0, ioj) yields e _1 (q — > 0, iu) = 1, rather than the correct e _1 (q — > 0, iui) = J°+ji , 
as shown in Fig. [3] where numerical results are compared with the Lindhard function. 

A similar discontinuity in GqWq corrections occurs for real metals such as Na and Al. 
Differently from the HEG, we find e _1 (0, iui) < 1, due to the nonzero interband contributions 



281 ] . In particular we obtain e Na (q = 0, 0) ~ 0.94, similar to the incorrect HEG value, and 



^(q = 0,0) ~ 0.008. The latter nears the proper Drude value, due to a substantial part 
of the aluminum Fermi surface being very close to a Brillouin-zone boundary, thus putting 
many of the metallic contributions of Eq. (fl3|) effectively into inter- rather than intra-band 
terms. For this reason, in the case of Al, the error induced by neglecting the intraband term 
is so small that the unphysical gap is almost invisible. 



IV. EXTRAPOLATED SMALL-q POL ARIZ ABILITY 

The solution we propose in this paper is devised to avoid the explicit (numerically ex- 
pensive) integration over the Fermi surface that would be required for a straightforward 
inclusion of the intraband term. We propose to compute the small q polarizability by a fit 
of the expected asymptotics. The time-reversal invariance implies the following symmetry 
of the matrix polarization: 

xoG,G'(q>^) = xo-c -c(-q,w) . (15) 

Therefore, the small-q expansion of Xoo,o(ci? includes only even powers. The expansion 
of the intraband term (J = j') in Eq. ffl3l) is 

Xoto(^^)^^-E ffkSfr-ejOc)) - q " Vk g |l + q-(MV k -.r|k,j)| 2 . (16) 
V B z ^ J ~ q • Vk£j 

The diagonal matrix element in Eq. (1161) is purely imaginary, therefore the last factor is 1 plus 
a q-quadratic contribution. For u = the intraband term is then a constant proportional 
to the density of states at the Fermi energy, plus quadratic corrections. For uj ^ the term 

8 



linear in q cancels because Vk^- is odd and the integral vanishes: the expansion begins with 
quadratic terms. The expansion of the interband j 7^ f term is easily seen to be no less 
than quadratic. To sum up, we use the following expression: 

xl t oM^) = ^ + Y,B^qrq s , (17) 

rs 

where A u , B^ s are real adjustable parameters, and A tuJp = for uj = iujp. The matrices B w 
are symmetric, and may have further symmetries depending on the crystal geometry. 

The off-diagonal elements G = G' 7^ of xo (the so-called "wings" of the matrix) are 
affected by a similar error, since they also contain the contributions of Eq. (JUJ). We also fit 
the intraband contribution to 

^W(^) = ^ G ' + E^ G '?r, (18) 

r 

where C wG ', D^ G> are complex adjustable parameters, and C KJpG = for u = iujp. 



V. RESULTS 



We determine the parameters A", B? s , C^ G ', and D? G ' in Eqs. (JUJ) and <$B) by a 
standard linear regression on values Xoo,o(q, and Xoo.c^ u ) computed for nonzero q- 
points inside a sphere of radius q c centered in T. We implement this procedure within the 



abinit 



12 



18] package. To test the effectiveness of the method, we apply it to the HEG 
in a simple cubic cell geometry, and to bulk sodium and aluminum in their experimental 
crystal structures (bcc a = 8.107 a , and fee a = 7.652 a , respectively). 

Figure H] displays the fitting of the computed polarizability xo of the HEG. Panels (a) 
and (b) compare the computed xo an d its small-q fitted parabolic expansion, for two differ- 
ent cut-off radii q c . Panels (c)-(f) display the resulting extrapolated small-q values of the 
polarizability as functions of the main parameters involved in the simulations and the fit. 
In these fits, the cut-off radius q c cannot be chosen too small, or else the number of q-points 
becomes insufficient to perform a reliable fit, especially at u) = 0, where the computed Xo is 
affected by significant numerical noise. Likewise, if q c is increased so much that it becomes 
comparable with the Fermi momentum kp, the outer points introduce a systematic error due 
to the non-parabolic q-dependency of Xo- An intermediate reasonably selected q c must then 
be adopted. This is slightly more important for u = iujp, where the fit is comparably more 
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Figure 4: (Color online) The parabolic polarizability Xo*) Eq. (fT71) . fitted to the computed 
Xooo(q, k>) °f the HEG (dots) restricted to q-points within in a sphere of radius q c centered at 
q = 0, and compared to the computed polarizability itself and to the exact (Lindhard) function, for 
(a) oj = and (b) oj = iojp. The computation involves a cut-off energy of 3 Ha, iVw = 16 x 16 x 16 
and a smearing temperature T smear = 0.005 Ha. Panels (c)-(f): convergence of the fitted values 
Xo t (qs ; w) (where the tiny q s = (7, 14,21) 10 -6 Oq 1 ) as a function of the cutoff radius q c , for differ- 
ent k-points sampling, and with (c) oj = and (d) oj = iojp, and for different smearing temperature, 
and with (e) oj = and (f) oj = iojp. Horizontal lines: the exact (Lindhard) values. 



sensitive to the value of q c , as shown in Fig H](d,f). Comparison to the Lindhard function 
show that this procedure provides a fairly accurate small-q Xo value, within few percent. 
Small smearing temperature is beneficial to a better accuracy in the determination of the 
asymptotic small-q behavior, but increase the numerical noise in the computed xo- 

As Fig. shows, the corrected screening successfully closes the unphysical gap. Of course 
in aluminum, where the fictitious gap is almost invisible, we see no significant difference in 
the GqWq corrections computed with and without the fit. The resulting curves are not very 
sensitive to the fit details, such as the value of q c , or A^ kpt . For the HEG we can compare the 
obtained bands with those computed via the Lindhard screening: the tiny almost uniform 
shift is due to the truncation in the number of empty states included in the summations of 
Eq. ( fl~3l) . which makes screening different in the large-q region, as illustrated in Fig. [3j 

Table |T] reports the occupied bandwidths of the metals studied in this work compared to 
previous calculations and experimental values. The comparison with the DFT-LDA values 
shows the well-known bandwidth reduction. The results for the HEG are close to Hedin's 
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Figure 5: (Color online) Comparison of the band energy obtained via the naive Go Wo calculation 
(dashed) to those obtained with the q ~ corrected polarizability (solid). For the HEG (a) the 
figure also shows Go Wo results obtained with the analytic (Lindhard) polarizability. For sodium 
(b) and aluminum (c) the bands are plotted along the (110) and (111) directions, respectively. 



0, 



computations [1], whi 
data in the literature 



e in the case of Na and Al the numerical value are comparable with 



29| | and approach the experimental values. 



VI. CONCLUSIONS 



In this work we have shown that the standard GqWq implementation of calculation of 
quasiparticle-corrected bandstructures, a basic tool to account for weak correlations in semi- 
conductors and insulators, describes metals correctly only in exceptional cases, like Al, where 
a substantial part of the Fermi surface falls very close to a Brillouin zone boundary, hence 
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Table I: Occupied bandwidth for the metals studied in this paper. The present results are compared 
to similar calculations and experimental data. For the HEG, the DFT-LDA result coincides with 
the free-electron model Fermi energy fp = En a (97r/4) 2//3 (r s /ao)~ 2 . All energies are in eV. 



r s /ao 


3 


HEG 
3.5 


4 


Na 
3.93 


Al 
2.07 


DFT-LDA 


5.57 


4.09 


3.13 


3.15 


11.01 


HEG G W a 


5.24 




2.86 






GqWq for metals b 








2.52 


10.0 


present work 


5.14 


3.84 


2.86 


2.81 


10.03 


experiment 








2.65 c 


10.6 d 



"Calculations by Hedin [1|J. 
b Calculations by Northrup et al. 29j|~ 
c Experiment by Lyo and Plumnieri [30 1 
^Experiment by Levinson et al. 31 1 . 



interband contributions make up for the missing intraband screening. In general (like in the 
HEG and Na), the incorrect intraband contribution to the small-q screening induces the 
opening of an unphysical gap at the Fermi energy. 

The proposed solution recovers the correct q — > polarizability by fitting a few small-q 
computed values, and solves this difficulty: the gap disappears, and the electron effective 
mass shows the expected few percent increase. This method requires a negligible computa- 
tional cost, contrary to other solutions based on Fermi-surface mapping. 

An entirely different solution can be devised, which avoids the limitations of the expansion 
(fl4"]) . and requires no fit altogether. Since the small-q polarizability needs to be computed at 
a q s much smaller than those generated by any practical k-sampling, it is possible to solve 
the KS equations on two k-point meshes shifted from one another by q s , and then apply 
directly Eq. ( Jl~3|) . We tried this method for Na and for the HEG, and find that the typical 
accuracies practically achievable in KS eigenvalues and wavefunctions allow us to compute 
Xo only for moderately large |q|, of the order of few percent of the k-mesh spacing. The use 
of such a not-so-small C[ s clS db representative of the q — > limit would induce systematic 
errors in the calculation of the GqWq corrections. The fit method is therefore practically 
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preferable. 



Acknowledgments 

The authors acknowledge S. Caravati for providing the routines for simulating the HEG 
within abinit. They also thank G. P. Brivio, R. Del Sole, P. Garcia-Gonzales, M. Gatti, 
R. W. Godby, V. Olevano, L. Reining, F. Sottile, and M. Verstraete for useful discus- 
sion. This work was funded in part by the EU's 6th Framework Programme through the 
NANOQUANTA Network of Excellence NMP4-CT-2004-500198. 



[1] L. Hedin, Phys. Rev. 139, A796 (1965). 

[2] L. Hedin and S. Lundqvist, in Solid State Physics, edited by F. Seitz, D. Turnbull, and 

H. Ehrenreich (Academic Press, New York, 1969), p. 1. 
[3] R. W. Godby, M. Schliiter, and L. J. Sham, Phys. Rev. B 37, 10159 (1988). 
[4] M. S. Hybertsen and S. G. Louie, Phys. Rev. B 34, 5390 (1986). 

[5] M. M. Rieger, L. Steinbeck, I. D. White, H. N. Rojas, and R. W. Godby, Comp Phys. Comm. 
117, 211 (1999). 

[6] G. Onida, L. Reining, R. W. Godby, R. Del Sole, and W. Andreoni, Phys. Rev. Lett. 75, 818 
(1995). 

[7] L. Reining, G. Onida, and R. W. Godby, Phys. Rev. B 56, R4301 (1997). 
[8] M. Rohlfing, P. Kriiger, and J. Pollmann, Phys. Rev. B 48, 17791 (1993). 
[9] M. Palummo, R. Del Sole, L. Reining, F. Bechstedt, and G. Cappellini, Solid State Comm. 
95, 393 (1995). 

[10] S. Massidda, A. Continenza, M. Posternak, and A. Baldereschi, Phys. Rev. B 55, 13494 (1997). 

[11] J. A. Soininen, J. J. Rehr, and E. L. Shirley, J. Phys.: Condens. Matter 15, 2573 (2003). 

[12] |http : / / www . abinit . org"7 

[13] http : //www. sax-project . org/. 

[14] |http: //www.f isica.uniroma2. it/~se lf /, 

[15] E. G. Maksimov, I. I. Mazin, S. N. Rashkeev, and Y. A. Uspensky, J. Phys. F: Met. Phys. 18, 
833 (1988). 



13 



[16] K.-H. Lee and K. J. Chang, Phys. Rev. B 49, 2362 (1994). 

[17] A. Marini, G. Onida, and R. Del Sole, Phys. Rev. B 64, 195125 (2001). 

[18] X. Gonze, G. M. Rignanese, M. Verstraete, J. Beuken, Y. Pouillon, R. Caracas, F. Jollet, 

M. Torrent, G. Zerah, M. Mikami, et al., Zeit. Kristall. 220, 558 (2005). 
[19] F. Bruneval, N. Vast, and L. Reining, Phys. Rev. B 74, 045102 (2006). 
[20] B. Holm and U. von Barth, Phys. Rev. B 57, 2108 (1998). 
[21] W.-D. Schone and A. G. Eguiluz, Phys. Rev. Lett. 81, 1662 (1998). 
[22] G. D. Mahan and B. E. Sernelius, Phys. Rev. Lett. 62, 2718 (1989). 
[23] E. L. Shirley, Phys. Rev. B 54, 7758 (1996). 
[24] S. L. Adler, Phys. Rev. 126, 413 (1962). 
[25] N. Wiser, Phys. Rev. 129, 62 (1963). 

[26] M. S. Hybertsen and S. G. Louie, Phys. Rev. B 35, 5585 (1987). 
[27] G. D. Mahan, Many Particle Physics (Plenum, New York, 1981). 
[28] Interband contributions are identically zero in the HEG. 

[29] J. E. Northrup, M. S. Hybertsen, and S. G. Louie, Phys. Rev. B 39, 8198 (1989). 

[30] I.-W. Lyo and E. W. Plummer, Phys. Rev. Lett. 60, 1558 (1988). 

[31] H. J. Levinson, F. Greuter, and E. W. Plummer, Phys. Rev. B 27, 727 (1983). 



14 



